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Abstract. Numerical evolution of the spherically symmetric, massive Klein- 
Gordon field is presented using a new adaptive mesh refinement (AMR) code 
with fourth order discretization in space and time, along with compactification 
in space. The system is non-interacting thus the initial disturbance is entirely 
radiated away. The main aim is to simulate its propagation until it vanishes 
near J?^ . By numerical investigations of the violation of the energy balance 
relations, the space-time boundaries of "well-behaving" regions are determined 
for different values of the AMR parameters. An important result is that mesh 
refinement maintains precision in the central region for longer time even if the 
mesh is only refined outside of this region. The speed of the algorithm was 
also tested, in case of 10 refinement levels the algorithm was two orders of 
magnitude faster than the extrapolated time of the corresponding unigrid run. 

Adaptive mesh refinement; numerical accuracy; numerical relativity 

1. Introduction 

The primary motivation to develop the presented techniques and perform the as- 
sociated investigations is the need for the efficient simulation of nonlinear dynamical 
systems. Numerical integration of nonlinear field equations is a difficult problem 
even if the metric is fixed. An obvious complication is that the field propagates in 
infinite space-time, but the computational resources are finite. Fortunately, infinity 
can be brought to finite distance by compactifying space-time, with a conformal 
rescaling of the metric pQ. A well-chosen coordinate transformation can also in- 
crease the efficiency of numerical calculations. However, a good transformation is 
often hard to find, especially when time evolution changes the system drastically, 
or in the presence of more than two different scales. Some examples are black hole 
formation, black hole merger, compact binary stars, etc. The sizes of such grav- 
itational radiation sources are very small compared to the produced wavelength, 
which is much smaller than the distance from the detector. Different length scales 
must be considered simultaneously, but it is extremely hard if time evolution is 
simulated on a uniform grid. Adaptive mesh refinement (AMR) algorithms over- 
come these difficulties by using a coarse base grid which is refined automatically at 
"interesting" locations for more precise calculation [2|. 

The precision also depends on the order of numerical schemes used. According 
to Hiibner, fourth order calculations are at least two orders of magnitude more 
efficient than second order Nevertheless, only second order calculations are 
known to be used by AMR codes in numerical relativity so far, although several 
implementations of the algorithm were developed since Choptuik's pioneering work 
0]. Our choice of the fourth order Runge-Kutta scheme is also supported by the 
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result of Hansen, Khokhlov and Novikov. They found that among the methods 
they investigated, this is the most efficient one in terms of accuracy and dissipation 
0. 

Problems where gravity is coupled to matters fields are complicated to start 
with, thus the code will be applied first to study simpler systems, physical fields in 
flat space-time. The simplest possible massive field is the free Klein-Gordon field, 
the time evolution of which is investigated in case of spherical symmetry in this 
paper. This system is known to provide certain surprises in numerical simulations 
UJ, moreover there is a straightforward means of checking the efficiency of the 
developed numerical method by making use of the analytic solution to the initial 
value problem. 

However, in contrast to the usual scenarios, there is no compact central object in 
this system that could keep any part of the initial disturbance from radiating away. 
Instead of studying the central region where vibration with decreasing amplitude 
is left behind, the aim is the simulation of the field at larger distances, where 
a hypothetical detector could be placed. As the wave propagates outwards and 
approaches its wavelength approaches zero, partly due to a physical effect [§], 
partly because of the space compactification. Hence for an accurate simulation of 
the long-range behavior, finer and finer subgrids are needed as the wave gets closer 
to J+. 

Interaction will be included in forthcoming studies. Plans include the verification 
of earlier numerical results on the logarithmic decay of <fi 4 breathers and the 
study of excitations of Bogomornyi-Prasad-Sommcrficld magnetic monopoles. In 
the latter case, a previous study [B] will be extended by the inclusion of the Higgs 
field's self interaction. In both of these settings, nonlinear massive fields evolve in 
fixed Minkowski space-time. 



2. The AMR algorithm 

The code is based on the Berger-Oliger algorithm .2 . At the initial time (T = 0), 
it starts with a uniform grid, where the values in the grid points are determined by 
the initial condition. Based on local error estimations of the first integration step, 
child grids may be generated and filled with values of the initial condition function. 
This procedure is applied recursively, until cither the local error becomes lower than 
the chosen tolerance in each point or the maximum refinement level is reached. In 
later integration steps, if a refinement level contains point(s) with relatively large 
errors, then finer levels are regridded in child to parent order (finest first). When a 
level is regridded, new points are interpolated using old points from the same level 
and the parent level. 

The same refinement ratio r is used on each level. For the refinement criteria, 
Richardson error estimation of an arbitrarily chosen function u is used. The grid is 
refined if 

\u(Q 2 AtA J{x,t))~u(Q 2AtaAx f{x 1 t))\ 
(1) ^TT— 2 > ' 

where QAt.Ax is the numerical integration operator, f(x, t + At) ~ QAt,Axf{x, t), 
q is the order of the method, / is the field and e is the error tolerance (see 0). The 
exact form of u proved to be irrelevant, provided that it has an approximate linear 
dependence on the field components. For example, the square root of the energy 
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r = 2, Runge-Kutta 4th order, n mar g in = 28 
t+2At 



t+At 



-28 -24 -21 -18 -15 -12 -9 -6 -3 12 

Figure 1. Subgrid margin narrowing in case of r = 2 and 4th 
order Runge-Kutta. Margin points -28, -1 at i, i + 2Ai, t + 4At 
etc. are calculated by interpolation of points in the coarse (parent) 
grid. 

density was a good choice in Klein-Gordon simulations, but the value of field / was 
equally good. For simplicity, I choose the latter in all simulations. 

To reach high precision in the numerical calculations, only symmetric finite differ- 
ence and interpolation schemes are used whenever it is possible. Space interpolation 
is simplified by aligning new subgrids with their parent grid. Time interpolation 
can be avoided in one space dimension, because it is possible to apply relatively 
large grid margins without noticeable efficiency loss. The initial subgrid margin is 
proportional to r and it shrinks in each step, until the time becomes equal to the 
time on the enclosing coarser grid. (See Fig.^) Then the new values on the margin 
are determined by space interpolation. 

Because of the above-mentioned proportionality, r should not be too large, oth- 
erwise the large margin sizes could result in unnecessary slowdown. In case of a 
larger refinement ratio, less refined levels are needed for the same precision, but a 
smaller r has the benefit that the mesh can adapt more closely to the solution. In 
my tests, simulations using double refinement were slightly more efficient (by about 
5%) than triple refinement, thus I choose r = 2. 

For time and space discretization, the same order of accuracy was used. Sec- 
ond order Runge-Kutta time integration with second order space discretization and 
third order space interpolation, fourth order Runge-Kutta with fourth order space 
discretization and fifth order interpolation. The presented numerical results are 
calculated in fourth order. A fourth order centered finite difference scheme uses 
5=2+1+2 points, thus the numerical error propagates 2 points in each step. How- 
ever, to avoid instabilities, an artificial dissipation term containing the sixth order 
sixth derivative was also applied y]: 

(2) af"""(Axf = a-(/ fc _ 3 -6/ fc - 2 + 15/ fe -i-20/ fe + 15/ fc+ i-6/ fc+2 + / fc+ 3)/A2:, 

where I used a = 0.01. It increases error propagation velocity to 3 points per step. 

Therefore refined grid margins must be shrunk by 3 points in each step. The 
fourth order Runge-Kutta method consists of 4 substeps, thus the overall loss of 
refined margin points is 12r = 24 in a coarse time step, see Fig. ^ However, a 
24-point wide margin would not be sufficient if the refined and the coarse grid have 
their origins (zero index points) or right edges at the same location. After a coarse 
time step, the coarse margin shrinks to 24-12=12 points. Although the leftmost 
refined margin point (with index i rc fincd = —24) coincides with the leftmost coarse 
margin point (i C oarse = —12) at this time, its neighbor (i ro fincd = —23) is not in 
a lucky position, it must be interpolated from coarse margin points. Fifth order 
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centered interpolation requires 3 points on both sides, but there is only one coarse 
point (icoaree = — 12) on the left side in this case. Choosing a slightly wider initial 
margin, with at least 28 points, can solve this problem. Then the coarse margin 
shrinks to 28-12=16 points which is enough to interpolate all refined margin points. 

For the fifth order interpolation of field values at point fc, the following formula 
is used: 

A Q f l0 + Aif tl + A 2 f l2 + A 3 f i3 + A 4 f u + A 5 f i5 



(•3) 



fk 



where io, «5 are the indexes of the six known points (from the coarse grid and 
the previous refined grid) and Aq, A§ are coefficients. The 6 coefficients are 
determined by fixing the function values, the derivatives and the second derivatives 
in points ii and i% . (Derivatives are approximated in fourth order of accuracy from 
2+1+2 points.) Their values depend on the distance (and existence) of the nearest 
old fine grid. The possible cases are shown in Table 

Table 1. Fifth order interpolation coefficients. 
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3. Simulation of a Klein-Gordon field 
A free, spherically symmetric Klein-Gordon field is described by the equation 

(4) d t 2 $ - -0£(r$) +m 2 $ = 0, 

where m is the mass parameter and Minkowski metric is used, ds 2 = dt 2 — dr 2 — 
r 2 da 2 . 

Fodor and Racz have shown that for arbitrary initial data with compact 
support the evolution of this system can be characterized by the formation of self 
similarly expanding shells built up of high frequency oscillations: 

(5) ®(t,r) w i" 3/2 * 4o (r/i) cos [m\/t 2 - r 2 + |) , 

where (f/t) are the contours of the shells; the form of this function depends 

on the initial condition. An important property of this solution is that both the 
frequency and the wave number of oscillations on the outer shell grows proportion- 
ally to y/t. To prove it, we determine these quantities as the partial derivatives of 
the cosine's argument: 

d{m\/t 2 — r 2 ) rat 



(G) 



dt yjt 2 - r 2 ' 



. . d(ra\/t 2 — r 2 ) mr 

The outer edge expands with the velocity of light, r w t, hence the frequency in a 
short distance Ar = t — r from the outer edge is approximately equal to the wave 
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number, their values can be estimated using the following formula: 

(8 us w k W ; w -^=. 

V(Ar)(i + r) \/2AF 

To remove the 1 /r coordinate singularity from Eq. 0] the "unphysical" field 
(9) /(t,r) = r*(t,r) 

is used instead of Space is compactified using a transformation as in Fiefs. [HI [0] : 



(10) T(t,r) = ut - y/oj 2 r 2 + 1, R(r) 



Vw 2 r 2 + 1 - 1 



wr 

where — oo < T < oo and < R < 1. 

To simplify numerical calculations, the second derivatives are removed from the 
field equation by introducing the partial derivatives of / as new variables fx and 
/r. Then the field equation (@J can be written as a system of 3 first order partial 
differential equations and a constraint condition: 

(11) d T f = f T , 

9tIt = - : , r, 9 /t — n/1 ^ n ,, -/r — %Rd R f T + — ^—dufR 



(12) 

(13) d T f R -- 

(14) f R - 

Boundary conditions. For the numerical calculation of derivatives near R — 0, 
information on the parity properties of the functions must be applied. Whenever 
the field $ is smooth, it has to be an even function of R. Hence / and fr are odd 
and fn is even: 

(15) f(-R) = -f(R), M-R) = -f T (R), fai-R) = fn(R). 

The behavior of function / near the other boundary is also important for similar 
reasons. The field $ vanishes in thus the function values for R > 1 points are 
supposed to vanish too: 

(16) f(R) = 0, f T (R) = 0, f R (R) = 0. 

The initial condition is a smooth, motionless hunch on the T(t, r) = hypersurface: 

if r > a — b and r < a + b 




(17) / 

(18) It 




otherwise, 



In the simulations, I used a hunch at R ~ 0.050 ± 0.037, with parameters u> = 0.05, 
a = 2 (center of hunch in r), b— 1.5 (half- width in r), c = 70 and d = 10. 

This disturbance is initially close to the origin, thus it is expected to reach J>^ , 
the future null infinity, at about T ~ 1 like the null geodesic denoted by dashed 
line on Fig. [21 
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f = arctan(0.015-(t+r)) - arctan(0.015-(t-r)) 



Figure 2. Penrose diagram of a part of the Minkowski space- 
time, with the i?, T coordinate lines. Initial condition is specified 
on the T = spacelike hypersurface (thick line). A null geodesic 
from the origin (R, T) = (0, 0) reaches at T = 1 (dashed line). 
A numerical solution preserves initial energy with 99% precision 
under its "validity limit" (upper thick lines). 



4. Results 

The simulations reproduced the most important feature of the analytical result 
by Fodor and Racz UJ, the self similarly expanding shells built up of higher fre- 
quency oscillations. The amplitudes of these high frequency "carrier waves" are 
modulated by the vf^ function in Eq. JSJ). As the shells expand and approach y + , 
the R-length (the wavelength in R coordinate units) of the carrier wave approaches 
zero, thus finer and finer grids are needed to simulate its propagation. Accord- 
ingly, to zoom into the vicinity of R = 1, the number of refinement levels must be 
increased, see Fig. 

At T = 0, the grid is only refined in a small range enclosing the initial hunch. 
Then the refinement follows the $ waves, it moves outwards and expands. When the 
waves approach the peak of the refinement level function begins to increase, 
see Fig. 0J Since the exact solution of the current problem is known, it might be 
possible to replace mesh refinement by a uniform mesh together with a clever time 
dependent coordinate transformation. The mesh refinement structure shown on 
Fig. 01 foreshadows the complications in finding such a transformation in space; a 
hunch like this cannot be smoothed using a simple, monotonous function for T < 1. 
There are complications also for T > 1, because the frequency of oscillations on the 
outer shell increases in time (see Eq. (JHJ). Moreover, even if an easily calculable 
transformation exists which makes mesh refinement unnecessary, the equal time 
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Figure 3. Exact (thin dashed lines) and numerically calculated 
field values (bold lines) on the T = 5 hypersurface. The top plot 
shows the result of a unigrid simulation with 256 grid points. The 
right edge is calculated more precisely, with 6 levels of refinement 
(middle). An even closer look with 10 levels of refinement is shown 
on the bottom plot. 




Figure 4. Number of mesh refinement levels as a function of R 
in different time slices. The top plot shows the refinement in T = 
0, 0.5, 1 (gray boxes), the bottom plot shows both the refinement 
in T = 5 (gray boxes) and the form of the $ function (black curve). 



surfaces would be different with the new coordinates, making it hard to transform 
the results back to the (T, R) system. 

Convergence. The 4th order convergence of the algorithm was verified by calcu- 
lating the time dependence of the convergence factor 
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Figure 5. Time dependence of the absolute error and the conver- 
gence factor in AMR and unigrid simulations. 



(19) Q ll^'-'-'ll 



||/Afl/2 - /II ' 

where 1 1 • 1 1 is the L 2 norm and /ah is the numerical solution of function / in case of 
a base grid with AR spacing. Locations and sizes of refined subgrids are stored and 
reused in calculations with the finer base grids (AR/2). Fig. [5] shows the absolute 
errors and the convergence factor for some unigrid and AMR simulations. The 
log 2 Q curves start with a plateau at a height of approximately 4, proving fourth 
order convergence. However, the order of convergence falls off near T ~ 1 because 
of the abrupt increase of absolute error when the wave reaches future null infinity 
(^ + ) and the i?-length of oscillations reaches zero, making even the finest mesh 
resolution insufficient. 

Energy conservation is violated numerically as the i?-length of waves approach- 
ing decrease below grid resolution. To check this, the following quantity is 
calculated instead of total energy: 

-Rmax T max 

(20) E(T max ,R max ) = J dRs(T max ,R)+ J dTj E (T,R max ), 

o o 
where e is the energy density and je is the energy flow. Numerically lost energy 
is mostly contained in the second term. By substituting -R max = 1, we get the 
total energy which is not conserved numerically. To "restore energy conservation" 
at a given T, i? max must be decreased. Both T max and i? max have a critical value 
below which energy is conserved with acceptable precision: E/Eq ~ 1, where Eq = 
E(0,1). Above the critical values, energy is lost: E/E Q < 1. Fig. shows the 
energy contour lines on the T max — i? max plane where only 95% or 99% of the initial 
energy is conserved: E/E = 0.95,0.99. The conservation bound can be pushed 
outwards both in i? max and T max by increasing the number of refinement levels. 
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FIGURE 6. Energy conservation bounds of numerical solutions cal- 
culated using different number of refinement levels and different 
error tolerance settings. 
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Figure 7. Time dependence of the error in the R < R ma , x = 0.375 
"central" range. 



Central range. As time evolves, the wave packet propagates outwards, leaving 
less and less matter in the center. The amplitude of oscillations decrease asymp- 
totically as i -3 / 2 oc T~ 3 / 2 (see Refs. [HI El)- Function / and its derivatives also 
decrease, thus one may suppose that a unigrid simulation is enough here. I tested 
this assumption by performing unigrid (256, 512 and 4096 points) and AMR simu- 
lations (256 points on base grid, 1 and 4 refinement levels). Errors and the norm of 
function / were calculated by restricting integrations to R < i? max = 0.375. Most 
of the matter leaves this range before T w f , thus function / becomes "smoother" 
and the grid is not refined at later times here. Consequently, the error of the 
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Table 2. Speed comparison of unigrid and AMR runs from T — 
to 5. The "n" columns contain the number of points on the 
(spacelike) base grid. Run times (to and ^amr) were measured on 
AMD Opteron (64 bit) 1.8 GHz hardware and Sun Java 1.5.0.02 
virtual machine. Nq and Namr are the total number of points 
of the spacetime grid. Tnstead of performing unigrid runs with 
32768 and more grid points, their times were extrapolated using 
the time of the n = 16384 run and assuming i un i gr id oc n 2 . 



unigrid 


AMR 












error tolerance 
















io- 8 


2 • 10~ 8 


io- 7 






n 


to 


n x 


levels 


*AMR' 


*AMR 


*AMR" 


N /N AM R 


to /t AMR 


256 


5.9s 


















512 


22.5s 


256 


X 


2 1 


18.8s 


17.6s 


17.3s 


1.7 


1.3 


1024 


73s 


256 


X 


2 2 


56s 


46s 


38s 


2.7 


1.6 


2048 


295s 


256 


X 


2 3 


118s 


108s 


84s 


4.5 


2.7 


4096 


1407s 


256 


X 


2 4 


293s 


263s 


194s 


7.7 


5.3 


8192 


5597s 


256 


X 


2 5 


630s 


586s 


468s 


12.9 


9.6 


16384 


22360s 


256 


X 


2 6 


1577s 


1416s 


1164s 


20.9 


15.8 


32768 


89440s* 


256 


X 


2 7 


4110s 


3574s 


2916s 


32 


25 


65536 


357760s* 


256 


X 


2 8 


10904s 


9590s 


8202s 


48 


37 


131072 


1431040s* 


256 


X 


2 9 


29064s 


26614s 21427s 


72 


54 


262144 


5724160s* 


256 


X 


2 io 


70148s 


63562s 47825s 


116 


90 



AR = 1/256 unigrid and the same base resolution AMR simulations are close 
for a while. However, the unigrid error starts to increase much faster at about 
T 12. This abrupt degradation proves that error can propagate inwards from 
outside (R > i? max ). Hence the unigrid error is only "acceptable" for T < 20, while 
the 4-level AMR calculation reaches the same level of inaccuracy much later, at 
about T w 150. See the curves with labels "256" and "256 x 2 4 " on Fig. □ On 
the same plot, it can also be seen that the error of a high precision unigrid run is 
smaller than that of the corresponding AMR with the same maximum precision. 
Compare the curve labeled with "512" (unigrid) to the "256 x 2 1 " curve (AMR, 1.5 
times faster), and "4096" to "256 x 2 4 " (AMR, 4 times faster). However, the AMR 
error curve increases much slower than the unigrid error curve, hence they meet at 
a certain time after which AMR is more accurate. If the AMR error tolerance is 
small enough, then the error at their meeting point is also small. 

Note that by restricting attention to the central range in i?, a much longer sim- 
ulation was possible in T. This feature can also be seen on the energy conservation 
curves on Fig. where T max increases when R max is decreased. 

Speed tests. The time of a unigrid simulation is proportional to the total number 
of grid points in the spacetime domain: 

(21) t run cx l/(Ax) d , 

where Ax cx 1/n is the grid spacing and d is the number of dimensions of the 
spacetime grid, d = 2 in this case. When using AMR, the same resolution can 
be reached much faster because only a small part of the grid is refined. In case 
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of 10 refinement levels, AMR would be two orders of magnitude faster than a 
corresponding unigrid run. AMR/unigrid run time ratios can be approximated by 
calculating the ratio of the total number of grid points in the spacetime domain 
(No for unigrid and JVamr for AMR), then adding the overhead (~ 30%) of the 
AMR algorithm. The £amrAo ~ 1-3-/Vamr/A?o formula is a good approximation for 
the measured times in Table |21 in case of 5 or more refinement levels. The formula 
(|21|l can be fitted nicely to AMR run times also. These fittings were performed at 
constant error tolerance values (10~ 8 , 2 • 10~ 8 and 10~ 7 ) for simplicity. The result 
is that for small T max values, the "effective dimension" d increases with T max until 
a plateau is reached near T max = 5, at a height of d = 1.35 ± 0.04. 

5. Other tests 

Further testing was performed with nonlinear 1 dimensional Klein-Gordon fields 
and periodic boundary conditions, $(x ± 1) = $(x), by adding a fourth order self 
interaction term to the Lagrangian: 

(22) L = 

The initial condition at t = is the motionless hunch used earlier in Eqs. (|17I18|) 
but substituting r by x, f by <fr, fx by $t and using the parameter values a = 0.5, 
b = 0.2, c = 10, d = 1. Test problems include massive (to 2 = 1, C4 = | and 
m 2 = 10, C4 = ^p), massless (to 2 = 0), free (04 = 0) cases and the sine-Gordon 
field: 

(23) L = i(9 t $) 2 -~(d x $) 2 +m 2 (cos$-l). 

In each of these cases, the initial hunch at x = 0.5 splits and its two parts propa- 
gate in opposite directions. At t = 0.5 they reach the periodic boundaries x = 0, 1 
where they pass through each other. They meet again at (t,x) = (1, 0.5), almost 
restoring the initial state but in a distorted form. As a consequence of the transla- 
tion invariance and the lack of nonlinear coordinate transformations, the derivatives 
and the numerical error do not increase substantially during the movement of the 
hunches, the number of required mesh refinement levels for a given error tolerance 
is determined by the initial state and preserved later. 

Another test was performed with the (p 4 breather described by the Lagrangian: 

(24) L(4>,d t <t>,d x <t>) = \(d t ^f - \(d x ^f -1(^-1)2. 

For the numerical simulation of this system, I compactified space using the same 
transformation as earlier in Eq. (|10f) . R = R(x) with uj = 0.05. Time t is not 
transformed. The initial condition is the same as that used by Geicke [7] , it contains 
a kink and an antikink in x = ±0.8. The radiation of the initial "kick" involves 
the outwards propagation of a wave with decreasing i?-length. As it propagates 
outwards, more and more levels of refinement are activated until it vanishes (in 
t > 10 3 if the base grid spacing is AR < 1/128). Then only the smooth, stable, 
oscillating wave remains near the origin and refinement is no longer needed. I 
performed simulations without further energy loss until t = 10 7 using a AR — 1/128 
base grid. However, I found energy loss in case of finer base grids. I also found that 
increased precision in the simulation of the initial radiation does not necessarily 
increase precision near the origin in the long run. The same amount of energy 
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disappears with the same rate if mesh refinement is not used at all in case of 
AR = 1/384 and t < 7000. Further experimenting is needed for the verification of 
the statement of Geicke [7] about the logarithmic decay of this system. 

6. Summary 

A new AMR code was developed for integrating field equations numerically in 
time. It is tested thoroughly using fourth order Runge-Kutta method and fourth 
order space discretization, but it is also possible to use other numerical schemes. 
The main test problem is the simulation of a spherically symmetric Klein-Gordon 
field in a special coordinate system l(TU)l with compactified space coordinate. The 
exact solution of this problem is known, thus it is possible to calculate the absolute 
error of the numerical simulation. By calculating the errors of AMR simulations 
with different base grid spacings, the fourth order convergence of the algorithm 
was shown. The numerical violation of energy conservation was also investigated; 
I determined the space-time boundaries of the "well-behaving" range. Inside this 
space-time volume, the sum of the total energy and its integrated outgoing flux is 
constant with an acceptable precision. The time boundary T max is a monotonically 
decreasing function of the maximum space coordinate -R max (see Fig.EJ. It means 
that energy is conserved numerically for longer time if a smaller central range is 
examined in space. The error of the simulation behaves similarly; longer runs can 
be "closer" to the exact solution if the error norm is calculated only in a small 
central range. The speed of the algorithm was also tested, in case of 10 refinement 
levels the AMR algorithm was two orders of magnitude faster than the extrapolated 
time of the corresponding unigrid run. 
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